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Abstract. We study analytically the equilibrium and near-equilibrium properties 
of a model of surfaces relaxing via linear surface diffusion and subject to a lattice 
potential. We employ the variational mean field formalism introduced by Saito for the 
study of the sine Gordon model. In equilibrium, our variational theory predicts a first 
order roughening transition between a flat low temperature phase and a rough high 
temperature phase with the properties of the linear molecular beam epitaxy equation. 
The study of a Gaussian approximation to the Langevin dynamics of the system 
indicates that the surface shows hysteresis when we continuously tune temperature. 
Out of equilibrium, this Langevin dynamics approach shows that the surface mobility 
can have different behaviours as a function of a driving flux. Some considerations are 
made regarding different dimensionalities and underlying lattices, and connections are 
drawn to related models or different approaches to the same model we study. 
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1. Introduction 

During the last decade, there have been great theoretical and experimental efforts to 
understand surface growth. This is due to possible applications, e.g., to the production 
of thin films and, from the basic point of view, to the interesting examples growing 
surfaces provide of non-equilibrium statistical systems IJ, in some cases with strong 
relation to relevant equilibrium systems 0. A very important example is provided by the 
discrete Gaussian (dG) model, which describes the universal features of the equilibrium 
roughening transition of many surfaces || . This transition is in the Kosterlitz-Thouless 
(KT) class, and thus the model is related to other important models featuring a 
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similar transition, such as the F model or the Coulomb gas |2|, p]. The dG model 
describes a surface minimizing surface area (to linear approximation), in which the 
surface height takes on integer values. Relaxing the latter condition leads to the sine- 
Gordon (sG) model for a real valued height field subject to surface tension and to a 
(lattice) potential favouring integer values of the field. The sG model is amenable to 
approximate analytic treatments 0, [5] which have allowed to develop a rather complete 
picture of the equilibrium roughening transition, and of the surface near-equilibrium 
properties as determined by Langevin dynamics || or kinetic Monte Carlo simulations 
0- 

There exist surface growth contexts, such as growth by molecular beam epitaxy 
(MBE), in which the most relevant relaxation mechanism taking place at the surface is 
surface diffusion Jl|. This, in turn, can be modelled as a surface minimizing curvature, 
instead of surface area. In reference ||, the following stochastic equation was proposed 
to study the interplay of this mechanism with a lattice potential favouring integer height 
values, similarly to the sG model 

~dt~ kAH ^ - ^f sm (-*r ) + " {r ' t] - (1) 

In ([!]), h(r,t) is the surface height above (a two-dimensional) substrate position 
r at time t; A is the two dimensional Laplacian and k, Vq, and a± are positive 
constants. F is also a constant representing, e.g., a driving flux of particles inducing 
the system to grow, rj is a Gaussian white noise with zero average and correlations 
(i](r,t)i](r' ,t')) = 2T5(r — r')5(t — t'), with T being temperature (we consider a unit 
Boltzmann's constant, ks = 1). In equilibrium [F = 0), this equation governs the 
thermal fluctuations of a surface described by (the continuum limit of) the Hamiltonian 

where i + 5 denotes a nearest neighbour to site i. In [0], numerical simulations of ([!]) 
showed an equilibrium roughening transition, similar to that in the sG model; namely, 
for temperatures below a critical value Tr the lattice potential is relevant and the 
surface is flat, whereas for temperatures higher than Tr the surface is rough. Out of 
equilibrium (F ^ 0), the surface mobility (to be defined below) behaves in different 
ways, depending on F and T. Although these simulations have been extended || 10 



no analytical approach had been made to study this model. In this paper we take a 
first step in this direction and apply to model (]l|)-(@) a variational mean-field approach 
successfully applied by Saito ]|] to the study of the sG model. In a subsequent paper we 
will refine this study by means of a dynamic renormalization group (RG) analysis of ([[]). 
The present paper is organized as follows. In section [| we briefly review the relationship 
between ©-(0) and related models for which approximate analytical and/or numerical 
results are available. In section [3] we study the equilibrium Hamiltonian (0) within the 
variational scheme of ||. Section [|is devoted to the approximate study of the Langevin 
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dynamics (JT[) within a Gaussian approximation for the probability distribution of the 
height. A discussion of the results obtained and our conclusions are found in section 
|5|. Some computational details on the solution of self-consistent equations relevant to 
section [3] can be found in appendix |6J], while appendix ^]2] discusses how the results 
are modified when considering the model on a triangular lattice (as opposed to the 
square lattice studied in the rest of the paper), and appendix |0] contains a discussion 
on results for substrate dimensions different from two. 



2. Background 



Crystal surfaces are often described within the solid-on-solid (SOS) approximation, in 
which the surface is characterized by a two-dimensional lattice height variable hi with 
hi/a± G Z and i being the square lattice position on a L x L dimensional substrate. 
Perhaps the simplest model is the discrete Gaussian model (dG) mentioned in the 
introduction, whose Hamiltonian is (we take a unit lattice constant) 

v 



n 



clG 



i+Sj 



hi/a± G Z 



(3) 



where v is a positive constant. Due to the difficulty to handle analytically the 
discrete sums in @, a continuum approximation is adopted introducing a potential 
that preserves the periodic symmetry in (|3]) and favours integer height values, leading 
to the sine-Gordon (sG) model (for reviews, see || |3|]) 

n sG = ^ h >- h ^ 2 + J2 v ° L f2nh 



1 — cos 



a± 



hi G 



(4) 



Both models undergo a KT-type roughening transition at a finite temperature Tjp 
between a flat and a rough phase. In the flat phase the roughness w 2 = (1/L 2 ) Ylifoi ~ 
h) 2 [where h = {^/ L 2 )(^2 i h i )} is finite and L- independent, while in the rough phase 
it diverges with the system size as w 2 ~ InL, i.e., as if we take Vq = in (^). In 
the flat phase, the lattice potential dominates and imposes a finite correlation length 
£. For T > TJp the lattice potential becomes irrelevant, although it modifies [H] 
the value of the surface tension u, and the correlation length diverges. Specifically, 
f ~ exp{C(T|, G - T)" 1 / 2 } for T -> Tf ~. 

As mentioned in the introduction, our aim is to study surfaces in which 
minimization of surface area is replaced by minimization of surface curvature, and we 
will thus replace model (ffl) by the Hamiltonian (0) proposed in ||. Interestingly, in the 



context of two-dimensional melting, Nelson [|TTJ] proposed (on the triangular lattice) the 
so called Laplacian roughening model 

2 



n 



LR 



E 



h; 



hi/a± G Z. 



(5) 



In equilibrium, numerical simulations and RG studies (see references in JT2|) indicate 
that © displays two phase transitions, both in the KT universality class, with an 
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hexatic phase between the two transition temperatures. However, the conclusions of 
analytical and numerical work by other authors (see a review in |Tj|) seem to be that 
model (|5]) features only one first order phase transition. Note model (0) is a natural 
continuum approximation of (|^) in the same spirit as the sG model is an approximation 
of the dG model. The numerical study || of Langevin dynamics ([[]) for (|2]) found an 
equilibrium continuous roughening transition between a flat phase and a rough phase 
in which w 2 ~ L 2 , the same behaviour of the so called linear MBE equation [that 
obtained by setting Vq = in (fj)]. Nevertheless, since the lattice potential modifies the 
value of the (in principle zero) surface tension, the long distance behaviour of the high 
temperature phase of model (0) was expected to be the same as that in the sG model. 
In what follows we apply Saito's variational treatment to @ and ([!]). In the case of the 
sG model, such mean field study @] allows to obtain the exact value of the roughening 
temperature, and a rather approximate estimation of the divergence of the correlation 
length near TJp. Thus we expect to obtain relevant information from this mean- field 
scheme. 



3. Variational mean-field method: equilibrium problem 

Following Saito [|J , our main assumption is that the most relevant features of model (fj) 
can be described by a simpler, solvable Hamiltonian: 

tto^X^te^M-*?) (6) 

where h(q) are the Fourier components of the height field 

3 

Here we consider periodic boundary conditions. Thus, q x = 2im x /L with n x = 
— (L — l)/2,...,L/2 and a similar relation holds for q y . Equation (^) defines a Gaussian 
Hamiltonian in which the values of S(q) are L 2 free parameters. We will fix them by 
minimization of the variational free energy Ty = !Fq+(TC— 7io)o> which is known to be an 
upper bound of the exact free energy T of model (0) by the Bogoliubov thermodynamic 
inequality [[J 

t<t v = T + (n-n )o (8) 

where JF is the free energy of model TLq and (• • -)o stands for the average with respect 
to the Boltzmann factor e~ n °^ T . 

Using the Hamiltonians @ and ([]) we obtain for the rhs of equation @ 



\ £ \n2nS(q) + \ Y^\q) ~ S~\q)]S(q) 

exp ( -r— w 2 ) } (9) 



2 q 2 

_ l_ exp f-- 



Variational theory of a model of crystalline tensionless surfaces 



5 



where we have defined S = T/{KUj{q)) with u(q) = 16[sin 2 (q , x /2) + sm 2 (q y /2)} 2 and 

" 2 = ^£"? = ^5>(9> ("» 

3 q/0 

(note that model (0) is symmetric under h — > — /i and thus, in equilibrium, ft, = 0). By 
minimizing jFy with respect to the parameters 5(g), we find they have to verify 



S~ (q) = .%Hq)+ l^-^e^f-^). (11 



7-1 /'^ _ C-1^„N , ,(_2 "0 _ „..2 



We can rewrite (|Tl|) by noting that the second term on the rhs does not depend on q. 
Hence 

s <«> = T^iw^ (12) 

where £ is a constant given by the self-consistent relation [note w 2 depends on £ through 
(0) and (HD] 

<" 4 = 4?r 2 V ° exp f-^w 2 ^) . (13) 
a 2 V 4 / 

Equations (|12|) and (|13|) are the solution to the equilibrium problem. We observe that 
the variational (Gaussian) approximation of Hamiltonian (fj) has a structure factor 
S(q) similar to that of the linear MBE equation. The only effect of the potential is 
to introduce a correlation length £ given self-consistently by equations ([13]) and fll2"l) . 
Among all mathematical solutions of equation fll3|) , the best approximation to model 
(@) is given by that value of £ that minimizes the variational free energy !Fy, which we 
denote by £ p h ys - Note all roots of equation (|T3|) can be easily shown to be stationary 
points of the function J-'viO- 

In order to proceed analytically, we need to take the continuum limit of the integrals 
appearing in @, fllQ|) . In this limit, we make the replacement L~ 2 J2 q ~^ (27r) ~ 2 J dq, 
and we can approximate co(q) = q 4 , hence using equation (|I2"D we get 
2 _ 1 /■ j _ T T£ 2 T 



w - J Aq ^ (q) + P) = t- ' w + °« »■ < 14) 

Keeping the dominant term in the above equation (in powers of £), and defining 
x = 2/t 1 / 2 a_i_T~ 1 / 2 7r~ 1 £~ 1 equation fllcf) becomes 

where 7 = QAVqo^kT" 2 ^ 2 . As shown in appendix |6.1| , there are different solutions of 
equation (|15|) depending on 7 (and therefore on temperature). Thus, = is always 
a solution of flTop, and is the unique solution for T > T c = 16Vq K^ 2 a±/&K. However, 
for T < Tq there appear two other finite solutions < £f x < £ 2 ^ 1 of equation fll5|) . 
In order to check which of the three roots provides £ p h y s in this temperature range, we 
compute the free energy difference 

^YLJ 1 s fT 2 ^ ~ Tv ^ 1 = 0)] " TeT " ^e-^/^D + 0(r 4 )- (16) 
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Figure 1. Variational mean free energy difference ATy as a function of the (inverse 
of the) correlation length for different temperatures. The values of £i and £2 are only 
displayed for the T < Tr case. The physical value of the correlation length £ p hys is 
given by the global minimum of AJ-y. For temperatures T > Tr the global minimum 
is always reached at = 0. Parameters used arc: Vq = a± = k = 1. 




We plot AjFy(^) in figure |T] for different values of T. For T < Tc, as can be seen 
in the figure, AjF^(^) has indeed a local maximum at ^jf 1 and a local minimum at 
£2 , while for T > T c both disappear. As derived in appendix for temperatures 
above T R = (e 1/2 /2)T c ~ 0.82 T c , the variational free energy difference has its global 
minimum at £ph ys — 0. However, for lower temperatures T < Tr, the finite correlation 
length £2 features a lower value of the variational free energy than the infinite correlation 
length solution, hence £ p h ys = £2 m this temperature range. Summarizing, within the 
variational approximation a roughening transition takes place at a temperature 

T R =^a ± ^V 1/2 - (17) 

Above Tr the correlation length is infinite and the surface is rough, with the same 
properties as the linear MBE model, i.e. S(q) ~ q A and w 2 ~ L 2 . Below Tr the surface 
is flat with a finite correlation length equal to £2- When we approach the roughening 
temperature from below, the correlation length does not diverge but, rather, tends to a 
constant value (see appendix |6.1| ) given by 



£(T - T 



Ana 2 N 1/2 



implying the roughening transition at Tr is of first order. Specifically, a cusp develops 
in the free energy Ty as a function of temperature at T = Tr, as depicted in figure 0. 

Although the results in this section have been obtained using a certain continuum 
approximation, we have numerically verified all our conclusions using the exact discrete 
sums in ([|) and (0). The exact variational results for the correlation length and the 
values of Tc and Tr for L = 1024 are compared in figure |2] to the approximate analytical 
expressions obtained in this section. We see that a first order transition indeed takes 
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Figure 2. Variational free energy Ty as a function of temperature for model (||) using 
the exact expression (^) (solid line) and our continuum approximation (dashed line). 
In both cases, Ty develops a cusp at T = Tr due to the jump in the physical value 
of £. Inset shows the values of Tr (lower curves) and Tc as functions of Vq within 
our continuum approximation (dashed lines) and using the exact discrete expressionis 

1/2 

(solid lines). In both cases, Tr cx Tc ~ V . Parameters used are a± = k = 1. 



place, although the values of Tr and Tc are modified. However, we still observe the 
non- linear dependence of Tr on Vo, see inset of figure |[ 

4. Dynamics within the Gaussian approximation 

In this section, we study the near-equilibrium dynamics of model @ by means of the 
generalized Langevin growth equation 

dhAt) 8H 

^r = F -JW) +v ' [t) (19) 

where T)i(t) is a white noise with correlations (T)i{t)r]j(~tf)) = 2T5ij5(t — t') and F is the 
flux of incoming particles in the surface growth picture, or a chemical potential difference 
in a generic context. This equation describes not only the non-equilibrium statistical 
dynamics of our model, but also the dynamics of the system fluctuations around the 
equilibrium state (for F = 0). Our approximation J| to the study of equation ([!]) will be 
to assume a Gaussian time-dependent probability distribution for the height field. Thus, 
we only have to calculate the first two moments of the probability distribution, namely 
the mean height^ h = (hi(t)) and the second moment (h(q,t)h(—q,t)) = S(q,t). Using 
equations @ and (^) (see a detailed account in [14]]) we find 

dh „ 2tcV q I . (2-Khi 



dt a± r (tt » } (20) 



dSM = -2TS(q,t) 



S (q)-S M + -^r{coB^—j 



dt 

% We make the homogeneity assumption that (hi(t)) is independent of substrate position, see 



(21) 
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where, within our Gaussian approximation, 

cos (M) \ = e-WW/.i cos (23) 

with u> 2 (t) being the time dependent surface roughnes. In all cases, we will study 
the set of coupled differential equations (p0|) and fl21|) subject to the initial condition 
/ij(t = 0) = for all substrate positions i. 

4-1. Equilibrium 

In equilibrium, i.e. for F — 0, the solution of equation (p0|) is ft, = [note (|22|) 1 and the 
solution of (pTD is the same as that of ([TT|) and (|13|) obtained in the previous section. 
The interest of equation (^j]) is that it allows us to study dynamically how does the 
system choose the physical value of the correlation length, and corroborate the results 
obtained in the previous section from the point of view of Langevin dynamics. Thus, we 
will integrate numerically the complete set of L 2 discrete equations ( f20| ) and (21) and 



perform the following experiment: starting from a flat surface and T = 0, we increase 
temperature by a certain (small) amount and wait until the system reaches equilibrium. 
Then, we increase temperature by the same amount and repeat the equilibration process. 
When the temperature is high enough (i.e. once the system is in the rough phase) we 
decrease temperature by the same amount and repeat the process of equilibration until 
T = is reached back closing a temperature cycle. 

We observe that the equilibrium first order transition found in the previous section 
indeed induces hysteresis in the system correlation length (see figure |3[) when the system 
is heated starting from T = 0, in the sense that the roughening transition takes place 
at the higher temperature T c and not at T R . The reason is that, for all T up to T c , 
the system stays in the local Ty minimum at £ 2 , even though for T R < T < T c the free 
energy already has its global minimum at = 0, since there is an energy barrier for 
the system to jump across the local maximum in AJ-'y Once the local minimum at £2 
disappears (i.e. for T > Tc), the surface is rough and exhibits an infinite correlation 
length. Conversely, when the system is cooled down starting at T > Tc, the system 
remains in the rough phase until T = is reached because £ _1 = is always a free 
energy minimum. 

4-2. Non- equilibrium 

In this section we allow F 7^ in ( p0|) and (fHJ) , in which case the former no longer has the 
trivial solution (h — 0). Rather, when the flux F is small (quasi-equilibrium condition) 
we expect the system to feature a structure factor S(q, t) of the same form as in 
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Figure 3. (Inverse of the) physical correlation length as a function of temperature, as 
determined from equations (|2|) and ( |l2| ) . The arrows indicate the heating and cooling 
experiment explained in the text. Parameters used are Vq = a± = n = 1 and L = 1024. 



equilibrium, all non-equilibrium effects reflecting in the (possibly non-trivial) behaviour 
of the average height. Actually, numerical simulations || [| of the full non- linear model 
(HI) seem to confirm this expectation. For this reason we neglect the feedback effect of 
the evolution of h(t) on the structure factor S(q, t) and take 

S(q,t)~ r YL_ 4 s (24) 

where £ is given by the physical equilibrium solution of section |||. Within this 
approximation, F C (T) = 2^±a exp{— 2tt 2 w 2 / a]_} becomes a constant and equation (^0|) 
can be written as 

dh „ „ 2nh , , 

— = F-F c sin (25) 

dt a± 

which is simple to integrate analytically (exact expressions for the solution can be found 
in [[|] and This equation has two different solutions depending on the values of F. 

If F < F c , then h tends to a constant value and the surface does not grow. If we define 
the surface mobility \i as 



1 /dh\ , , 

" = (26) 

where the overline stands for average over a time larger than /i -1 , then for F > F c one 
obtains from the exact solution of (p5|) a non-zero value for fi: 

F _2N 1/2 



"=^-?IJ ' ( 27 ' 

In figure ^| we plot the surface mobility as a function of T. Using the equilibrium 
solution for £ described in section |3], for temperatures above roughening (T > Tr), we 
have that £~ h y S = 0, which implies F c = and /i = 1. Thus, above roughening the surface 
shows linear growth with a maximum (unit) mobility. In the flat phase (T < Tr) the 
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Figure 4. Surface mobility as a function of T for different values of the driving flux 
F. The values of the mobility are obtained from equation (26) using equations ( |l2|) 
and ( |l3| ) with the parameter values Vq = k = a± = 1. 



mobility is equal to zero (i.e. the surface does not move) for a small flux F < F C (T). 
For larger values of the flux (F > F C (T)), the mobility depends nonlinearly on T for 
all temperatures up to T R . Due to the jump of the correlation length at T = T R , 
the mobility also has a jump at this temperature value. These three behaviours of the 
surface mobility as a function of temperature and driving flux agree with those obtained 
H H for the full model (|l|), except for the discrete jump of \i at T = T R . 



5. Discussion and conclusions 



Summarizing the equilibrium results obtained in the previous sections, the variational 
approximation predicts a first order phase transition for model (0), and the associated 
hysteresis phenomenon. In particular, within a Gaussian approximation, Langevin 
dynamics predicts that a rough surface can preserve its infinite correlation length when 
cooled down across the roughening temperature. Moreover we have found that these 
results apply both on the square and on the triangular lattices (see appendix |6.2| ). 
Hysteresis behaviour and a first order transition have been reported in |H| and [|l(J 
for models related with the Laplacian roughening (LR) model @ on the square lattice 
and for the LR model on the triangular lattice |17| . However, as mentioned in the 
introduction, other authors seem to obtain two KT transitions for the LR model both 
on the square [31| and on the triangular [|l9j lattices. Note that our Langevin dynamics 
results within the Gaussian approximation yield a discrete jump in the surface mobility 
\i at Tr, which is not found in simulations of the full nonlinear model (|]). This might 
indicate that the first order character of the transition is in our case an artifact of 



the variational approximation. Moreover, this approximation (see appendix |6.3|) also 
predicts a phase transition for model (0) in d — 1, which is also obtained for the 
sG model. This result points out the limitations of this approximate framework for 
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situations in which fluctuations are very relevant for the system behaviour (as in the 
d = 1 case). Since model (|2|) features intrinsically strong fluctuations (as does e.g. 
the linear MBE equation [[]), it is desirable to go beyond our present mean field 
approach to this model. We can take two steps in this direction. One (numerical) is 
to perform extended simulations of both the LR model and model (|2|) [or, equivalently, 
its equilibrium Langevin dynamics (ffl)]. The results \W\ seem to indicate that in both 



cases there is only one continuum transition, though with strong size dependence in the 
LR case for sizes up to moderate (but not large). The other (analytical) improvement 
is to perform a dynamic RG analysis of ([!]) along the lines of f| for the sG model. This 
study is particularly important bearing in mind that the lattice potential is expected 
to contribute a surface tension, absent in equation ([!]), which should then dominate 
the scaling behaviour as compared with surface diffusion Jl|. This phenomenon is 
clearly beyond our mean-field approach, which neglects parameter renormalization. 
Moreover, comparing with the sG case, in the latter the variational mean field [[| and 
the perturbative (in powers of Vb) RG || approaches yield the same (exact) roughening 
temperature T|P = 2va\/ir, which is independent of Vq. However, in our case Tr does 
depend on the lattice potential (and on k) as a fractional power. Thus, we do not expect 
a perturbative RG treatment to quantitatively agree with the expression for Tr derived 
here. This will be the subject of a forthcoming publication. 
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6. Appendix 

6.1. Solution of the self- consistent equations 

In this Appendix we calculate the self-consistent solution of equations ([12]) and ( p~3| ) for 
the equilibrium correlation length of the variational approximation (^|) to model (|2]). 
By defining x = l^a^T^I 2 ^^ 1 and 7 = 64^1 kT"^- 2 , equations (|T2|) and (|IJ) 
become, within the continuum approximation made in section § 

It is obvious that equation fl2"8[) always has the solution x = 0, and that for some values 
of 7 it may also have non-zero solutions. Our first aim is to determine the critical value 
of 7 for which x = is the unique solution. To this end, we rewrite the equation in the 
following way 

x = 7 VVi/4< (29) 

Now the solutions are the intersections of the function y = f{x) = 7 1 / 4 e - 1 / 4;c2 with the 
straight line y = x. As we can see in figure [|, for 7 > 7c there are three solutions of 
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Figure 5. Graphical representation of equation (J29|). The dashed line is the y = x 
function, while the solid lines show y = -y 1 / 4 e _1 / 4:E for different values of 7. 



equation (^9|), two solutions for 7 = 7c and only the trivial solution x = for 7 < jc- 
The value of 7c can be calculated using that for 7 = 7c the unique solution x = x s 7^ 
verifies (^) and also the equation 

1 



1/4 



(30) 



obtained by requiring that the slopes of y = x and y = f(x) be equal at x = x s . With 
these two equations it is easy to obtain x s = 2~ 1//2 and 7c = e 2 /4. Using the definition 
of 7, the temperature for which = is the only solution of equation (|T3| ) is then 
given by 

16V 1/ \ 1 / 2 a ± 



■ c 



C7T 



(31) 



Now, for T < Tc we have to determine which of the three solutions of (^8|) provides the 
physical correlation length. Since £ -1 = is the unique solution for high temperatures, 
we take as a reference value JFy(^~ 1 = 0), and note that AjF(^) = jFy(^)— jFy(^ _1 = 0) is 
stationary at any root of equation ([13]). Thus, we will consider as the physical solution for 
the correlation length that root of (^) for which Ty has an absolute minimum. Starting 
out with high temperatures, the condition AjFy(^) = will signal the temperature at 
(and below) which = ceases to be the global minimum of the variational free 
energy and thus the system physical correlation length. Using our previous notation, 
the condition AjFy(^) = reads 

x 2 = y e -V- 2 (32) 

where 7' = 64/?a 2 L Vo/(7r 2 T 2 ). Using the same argument as above, it is easy to show that 



for 7' < 7^, = e there are non-zero solutions of fl32"l) . This means, using the definition of 
7', that there is a temperature given by 

e l/2 

T R = —T c , (33) 
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such that for T <Tr the global minimum of the free energy is attained for a correlation 
length £ 7^ 0, whereas for T >T R the physical solution is £~ h y S = 0. 

6.2. Triangular lattice 

The Laplacian roughening model was initially proposed by Nelson on the triangular 
lattice |[T|. Thus, it is worth studying how do the features of our model (0) change 
when the substrate geometry is different from the square lattice considered in the 
text. Nevertheless, we expect that only nonuniversal quantities — such as the transition 
temperature and the numerical value of correlation length — depend upon the lattice 
geometry. The Laplacian roughening model on the triangular lattice is given by 

2 



/Xhj - h i+s ) 



8 



(34) 



with i + 5 being any of the six nearest neighbours of site i. For this case [14 
co(q) = lQ{sin 2 (q x /2) + sm 2 [(q x + y/3q y )/4] + sm 2 [(q x - ^/3q y )/A}} 2 , where q = 
(n x /L)b x + (n y /L)h y , with = 2n{e x — (l/v3)e 2/ ] and h y = (An/ y/3)e y , where 
Hi = — (L — l)/2, ...,L/2 and ej are the standard basis vectors. In the continuum 
limit, So (9) — 4T/(9k<7 4 ), and we recover equation (0). Taking the continuum limit 

(i-e E q -*• 4 Jbz (£f - *F Jo (2/v/5)1/V ^ d ^' where BZ denotes the first Brillouin 
zone), we get 



2 T 

w ~ 



2ttk(2/V3) 2 Jo g' 4 



-dq'. (35) 



(2/V3) 2 

Thus, by defining T' = T/3, and £' = 3 1//4 £, we get the same equation (|T3D but with 
redefined constants T' and One can readily reproduce all the results obtained in the 
text, simply by making the replacements T — > T' and In conclusion, on the 

triangular lattice a first order roughening transition is also obtained, the only effect of the 
geometry being a shift in the value of the roughening temperature 7 1 ^ riang - = T^ quare /3. 

6. 3. Substrate dimensions d 7^ 2 

In this appendix we discuss the possibility to find a roughening transition in equilibrium 
when model @ is defined on a substrate of generic dimension d. In that case, equation 
(O) is still valid, but with 



^ g / fe<PFi ' (36) 

within the continuum limit. For substrate dimension d > 4, the integral (^) is finite 
for £ -1 = 0, namely w 2 ^ 1 = 0) = fQ7r d_4 T '/[/c(<i — 4)] (where -fQ is the ci-dimensional 
angular integral K d = J d d - 1 Q/(2n) d = 27r d / 2 /[(27r) d r(d/2)]). Thus, = is no 
longer a solution of equation (|T3|). Therefore the system has no rough solution and is 
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in the flat phase for all temperatures. On the other hand, for d < 4, the integral above 
may be approximated by 



^ 2 (0 < - , ■ ~r~4 /a\ e 4 - d . (37) 

k 4 sm(7ra/4j 

In this case, = is always solution of equation ([jjD, there being two additional finite 
solutions when T < Tq. The value of can be calculated using the same argument as 
in d = 2 and is 

d<4 _ 8nal sin(d7r/4) / 4tt^\ ^ 
Jc " (4-d)K d e7r3 { K al J ' ^ 

In order to know which solution of equation (|T3|) minimizes the variational free energy, 
we calculate AjFy, which now reads 

TL rf d 8sin(7rd/4)^ T 1 ' 

In this case, the local minimum ^ is also the global minimum and the physical 
solution for temperatures below the roughening temperature (T < T^ <4 ), which is given 
by 

T^ = ^T d c <\ (40) 

For temperatures above roughening (T > T^ <A ), = provides the global free energy 
minimum. Thus, for d < 4 there is a first order roughening transition at T^ <4 . Note this 
includes d — 1, which might seem conflictive since in this case model (|2]) is expected 
to be in the rough phase for all values of T (|TTJ: In d — 1 thermal fluctuations are 
expected to destroy the ordered flat phase for any temperature value. Our result can 
be understood by noting that Ty is not a true free energy, in the sense that it is not the 
free energy of any modelQ, but rather an upper bound for the free energy of model @. 
Actually, one obtains exactly the same result in the variational study of the sG model 
in d = 1 H . Note that in this reference the analysis of the d = 1 case is incomplete, 
with the incorrect conclusion that the variational theory predicts no phase transition 
when d = 1. The complete expression for AjFy(^) analogous to (|39|) indeed shows that 
also for the sG model in d = 1 the variational approximation does predict a non-zero 
temperature below which the physical value of the correlation length is finite. 

Finally, for d = 4 equation ([13]) is very similar to that obtained by Saito for the 
sine- Gordon model 

< -4 = 4 l!^( 1 + 7r 4 e 4 ) -Ti£ T> (41) 
CL 1 

Following Saito's analysis for the sine-Gordon model 0, we readily obtain that for 
d = 4 our model has a Kosterlitz-Thouless transition when T = T^ =4 = \Qna\. The 
correlation length now diverges as £ ~ exp{— A/(T — T^ =4 )} when T — > T^ -4 ' - (A is a 
T independent constant). 

+ The free energy of model Ho is Tq. 
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In summary, within the variational approach, our model displays a first order 
transition for d < 4 between a flat phase and a rough phase with the properties of 
the linear MBE equation. For the marginal dimension d = 4 this transition becomes of 
the Kosterlitz-Thouless type, whereas for d > 4 the surface is in the flat phase for all 
temperature values. 
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